Penetration and Displacement Behavior of N2 in Porous Interlayer Structures Containing Water/Salt Component by Molecular Dynamics Simulation

The penetration and displacement behavior of N2 molecules in porous interlayer structures containing a water/salt component with porosities of 4.29%, 4.73%, 5.17%, 7.22%, and 11.38% were explored using molecular dynamics simulations. The results demonstrated that the large porosity of the interlayer structures effectively enhanced the permeation and diffusion characteristics of N2. The water and salt in the interlayer structures were displaced during the injection of N2 in the porosity sequence of 4.29% < 4.73% < 5.17% < 7.22% < 11.38%. The high permeance of 7.12 × 10−6 indicated that the interlayer structures with a porosity of 11.38% have better movability. The strong interaction of approximately 15 kcal/mol between N2 and H2O had a positive effect on the diffusion of N2 and the displacement of H2O before it reached a stable equilibrium state. The distribution of N2 in porous interlayer structures and the relationship between the logarithm of permeability and breakthrough pressure were presented. This work highlighted the effects of porosity on the permeability and diffusion of N2/H2O in the interlayer, thus providing theoretical guidance for the development of petroleum resources.


Introduction
Oil has aroused growing attention as a crucial energy resource around the world [1][2][3]. Oilfield exploration technology is vigorously investigated to discover more oil resources [4][5][6]. The porous interlayer is one of the most important factors affecting the distribution of the remaining oil, which increases the heterogeneity of the reservoir [7], complicates the oil-water movement relationship [8], and disorders the distribution regularity of the remaining oil [9,10]. Previous studies have mainly focused on explaining the influence of the interlayer on oil-gas movement [11]. By exploring the effect of three types of interlayers on oil-water movement and the formation of rised water, Wang et al. [12] found that the oil layer can be divided into several connected or disconnected units and that it controls the movement of oil and water and ultimately affects the distribution of remaining oil. Luo et al. [13] discovered that the vertical migration interval of the oil and gas in the heterogeneous transport layer is hindered by the interlayer, thus complicating the migration path. Wu et al. [14] found that the interlayer restricts the vertical permeability of the reservoir. These results indicated that the interlayer has an important influence on oilfield exploration technology [15]. Considering the complexity of experimental conditions, theoretical calculation provides an effective way to explore the interlayer. The process of gases, including CO 2 and N 2 , passing through the interlayer is the best method to define structural features. Fang et al. found that slug injection of CO 2 was a promising technology for the enhancement of oil recovery in an MD simulation. However, few studies have focused on the injection of N 2 into the interlayer [3], which could effectively open the channel for the migration of oil and gas and cause the remaining water/salt component to be the only substances in the interlayer. Therefore, it is vital to elucidate the effects of porosity on the interlayer for N 2 penetration, H 2 O displacement behavior, and the diffusion mechanism.
In this work, we investigated the penetration and displacement behavior of N 2 in interlayers containing a water/salt component and with porosities of 4.29%, 4.73%, 5.17%, 7.22%, and 11.38% using molecular dynamics (MD) simulations. First, porous SiO 2 was built to simulate the interlayers. Subsequently, the pores of SiO 2 were filled with a salt solution, followed by the displacement of N 2 molecules under breakthrough pressure. Next, the diffusion process and mean squared displacement were calculated in consideration of the porosity and pressure effects. Then, the gas distribution and the permeability-pressure relationship were analyzed to show the microscopic mechanism of N 2 /water/salt passing through the interlayer. Our results highlight the potential impact of the interlayer on the gas diffusion process, providing theoretical guidance for oil exploitation.

Model and Method
The silica structure was selected as a typical mineral component [16]. The generation of the frameworks was conducted using a 10 × 10 × 10 supercell of silica structure with dimensions (xyz) of 71.60 Å × 71.60 Å × 71.60 Å. Considering that the remaining oil exists in interlayer structures with porosities of 4-8% [17], pores along the z-axis with porosities of 4.29%, 4.73%, 5.27%, 7.22%, and 11.38% were established in the center of the silica structures. These structures were tagged as P4.29%, P4.73%, P5.27%, P7.22%, and P11.38%. Furthermore, the silica structures used in this study were rigid and fixed, and the silica surface of the pore was hydroxylated and set at 9.6 nm −2 [17].
MD regards all matter as a particle system composed of atoms and molecules, and the motion of all particles complies with the laws of classical mechanics or quantum mechanics. The most important characteristic of MD is that it ignores the interactions of atoms between electrons and electrons and the electrons and the nucleus, only regarding the atom as a basic unit to calculate the energy of the whole system and thus greatly improving the computational efficiency, making the MD method applicable to larger computing systems. Therefore, the MD simulation was adopted to complete these calculations. A saline solution with an H 2 O to salt molar ratio of 1000:1 [18] was uniformly distributed in the hydroxylated silica pores to simulate the formation of brine. An N 2 molecule box was placed on the left side of the silica structures in order to create a pressure differential in the initial model. A vacuum slab was placed on the right side of the silica structures. A rigid graphene sheet was placed on the left side of the N 2 flooding and the right side of the vacuum slab as a rigid piston (Figure 1). For the system, the pores were uniformly filled with water/salt at the initial moment. A 1 ns equilibrium MD simulation with the NVT ensemble was conducted to achieve an equilibrium state [19]. The well-known 12-6 Lennard-Jones potential 4ε[(σ/r) 12 − (σ/r) 6 ] [20] was used to calculate the van der Waals interaction, where r represents the distance between two atoms. The cutoff distance was set to 11 Å. The rigid SPC/E model [21] was used for modeling water molecules with the SHAKE algorithm [22], while the Lennard-Jones potential was adopted to model the neutral N 2 [23]. The silica surface was modeled by the CLAYFF [24] force field. All MD simulations were performed in the NVT ensemble (the number of molecules, volume, and temperature were constant) with a time step of 1 fs [25]. The temperature, kept at 298 K, was controlled by a Nosé-Hoover thermostat [26], and every 1000 steps, a frame was set for data analysis. LAMMPS software [27] was employed, and total simulation was conducted for 1 ns. The application of the periodic boundary condition was done in all directions [28]. The desired transmembrane pressure in the experiments was applied to the feed piston to push the water/salt across the membrane.

Discussion
During the injection process, N 2 passed through the pores of SiO 2 with breakthrough pressures of 35.83, 14.00, 11.18, 5.77, and 0.50 MPa, which were selected to put the interlayer porosities of 4.29%, 4.73%, 5.17%, 7.22%, and 11.38% in the displacement systems, according to the experimental data (offered by the Research Institute of Exploration and Development, PetroChina Xinjiang Oilfield Company). Breakthrough pressures were applied to the rigid graphene sheet. Figure 2 presents the process snapshots of N 2 injection in the porosities of 4.29%, 4.73%, 5.17%, 7.22%, and 11.38%. For a certain porosity, H 2 O, Na + , and Cl − ions were displaced after the injection of N 2 , and the displacement degree of H 2 O, Na + , and Cl − gradually increased with time. Taking the structure with a porosity of 4.29% as an example, all N 2 molecules were gathered on the left side of the SiO 2 initially. After 0.3 ns, some N 2 molecules passed through the SiO 2 phase, and meanwhile, some H 2 O, Na + , and Cl − were displaced into the vacuum phase and the N 2 phase. With an increase in simulation time, more N 2 , H 2 O, Na + , and Cl − passed through the SiO 2 phase until the system reached equilibrium at 0.7 ns. In comparison, the larger porosities led to more evident displacement phenomena after a period of time. When the porosity was 11.38%, almost all the H 2 O dissolved into the N 2 phase. Detailed observations indicated that after injecting N 2 , all components were in dynamic equilibrium within the interlayer after 1 ns. To analyze the displacement process, the center of mass (COM) of N 2 was explored, as shown in Figure 3. For certain interlayer structures, the COM of N 2 moved toward the right of the interlayer structures. After 0.3 ns, the COM of N 2 moved −0.65 Å, −2.52 Å, −1.30 Å, −2.24 Å, and −1.93 Å for P4.29%, P4.73%, P5.27%, P7.22%, and P11.38%, respectively. In other words, the COM of N 2 moved far toward the left of the interlayer structures at the beginning of the simulation, except for the structure with a porosity of 4.29%. This was attributed to the achievement of an equilibrium state. In the next 0.7 ns, the COM of N 2 moved toward the right of the interlayer structures with moving distances of 5.66 Å, 5.87 Å, 7.74 Å, 7.14 Å, and 16.29 Å for P4.29%, P4.73%, P5.27%, P7.22%, and P11.38%, respectively. After the same simulation time, large porosity led to an increase in the movement distance of the COM, which was consistent with the process snapshot analysis. For P11.38%, the large pore size led to more gas molecules passing through interlayer structures, resulting in a faster passing process than other structures. Therefore, the curve of P11.38% shows an extremum and a large slope.

Movability Evaluation
Permeance is an indispensable indicator to determine the displacement efficiency of N 2 over water/salt [3]. Table 1 shows the permeances of N 2 passing through SiO 2 with five interlayer structures with different porosities. As the porosity increased, the permeance of N 2 gradually increased in the sequence of P4.29% < P4.73% < P5.27% < P7.22% < P11.38%. A large pore size results in a small barrier, resulting in a small N 2 energy repulsion. Moreover, the weak interaction guarantees that the gas molecules pass through the pores of SiO 2 rapidly. In order to evaluate the diffusion behavior and movability of N 2 and H 2 O molecules, mean squared displacement (MSD) [29,30] was calculated: where R i (t) represents the position of i atom at t time, and R i (0) represents the initial position. Each pore was divided into eleven layered blocks, and MSD was calculated for each structure. Previous studies confirmed that the gas diffusion rate was positively correlated with the slope of the MSD profile. As exhibited in Figure 4a, on the whole, the MSD of H 2 O along the z-axis in the five interlayers increased sharply at the beginning of the simulation (0-200 ps). The slope of the MSD followed the sequence of P11.38% > P7.22% > P5.27% > P4.73% > P4.29%. In particular, the structure with a porosity of 11.38% showed the fastest gas diffusion because the large pore size enhanced the H 2 O diffusion rate. After the initial increase in the MSD, the slope of the MSD flattened out to a constant value because of the equilibrium state of the system. During this period (200-1000 ps), the equilibrated MSD of H 2 O followed the sequence of P11.38% > P7.22% > P5.27% > P4.73% > P4.29%. The trend of the MSD curve agreed well with the physical characteristics and porosity of the pores.  During the simulation, the porosity of the interlayer structures increased from P4.29% to P4.73%, P5.27%, P7.22%, and P11.38%, and the equilibrated MSD of N 2 increased to 0.17%, 33.44%, 45.11%, and 117.61%. This result confirms that the MSD of N 2 is more significantly enhanced by the large porosity of P11.38% than by those of other structures. Conversely, the MSD of H 2 O is more affected by porosity. In short, increasing porosity can enhance the MSD of H 2 O more than that of N 2 at both early and late stages of the simulation, particularly for the large P11.38%. Therefore, it can be envisioned that large porosity can greatly promote the diffusion of N 2 in SiO 2 for the displacement of water and salt. Higher mobility of the gas molecules in the middle of the interlayer was observed for structures with large porosity, and the moving ability near the small silica surface was poor. However, over time, the overall effect of the H 2 O slug flooding became more severe than that of N 2 .

Distribution of H 2 O, N 2 , and Salt in Interlayer Structure
In order to evaluate the distribution of the H 2 O, N 2 , and salt in the Z direction visually, the relative concentration profiles over 1 ns were calculated, as shown in Figure 5. On the z-axis, the range of 0-71 Å represents the N 2 phase, the range of 75-165 Å represents the porous SiO 2 phase, and the range of 170-240 Å is the vacuum phase. The rigid graphene sheet was placed on the left side of the N 2 flooding and the right side of the vacuum slab as a rigid piston. Initially, most of the H 2 O, Na + , and Cl − was in the pores of the SiO 2 . After 100 ps, the density difference of the N 2 /H 2 O/salt at the SiO 2 interface gradually diminished, while a substantial change in the density of the N 2 /H 2 O/salt inside and outside the pores of the SiO 2 took place. The small pore size of the SiO 2 with a porosity of 4.29% limited the number of molecules in the free phase, resulting in almost no N 2 molecules left in the pores. The lower porosity led to a smaller volume to fill. And less N 2 in the free phase. However, some H 2 O with dipole moment, Na + , and Cl − ions existed in the pores of the SiO 2 , which is explained by the larger porosity. On the whole, some N 2 molecules can pass through the pores of the SiO 2 into the vacuum phase in all five systems under the breakthrough pressure. However, most N 2 molecules remain in the pristine N 2 phase. The high peak in the range of 0-10 Å was attributed to N 2 adsorption in the graphene sheet surface, and another peak in the range of 75-80 Å was attributed to N 2 adsorption in the SiO 2 surface, while the remaining N 2 molecules were uniformly distributed. The H 2 O and salt diffused from the pores of the SiO 2 into the N 2 and vacuum phase, and the number of H 2 O, Na + , and Cl − ions in the N 2 phase was smaller than that in the vacuum phase. The results show that N 2 slugs can effectively displace the saline solution in the pores of the SiO 2 into the vacuum phase.

The Relationship between the Parameters
To investigate the relationship between permeability, porosity, and breakthrough pressure, the proper function fittings are discussed. Figure 7a shows the relationship between the permeability of N 2 and the breakthrough pressure of N 2 . After logarithm fetching on the permeability of N 2 (Log(K)), a linear relationship with the logarithm fetch breakthrough pressure of N 2 (Log(P)) was given. The fitting formula is Log(P) = −0.775 × Log(K) − 3.298, which is in agreement with the experimental results of Log(P) = −0.658 × Log(K) − 0.557, which was called by log-log analysis. The same order of magnitude for the line slope and the experimental results shows the rationality of the model and accuracy of the calculation parameters. Furthermore, the large permeability can be obtained by adjusting the breakthrough pressure according to the fitting formula. In addition, the intersection profile of the porosity and gas injection breakout pressure was fitted, indicating an exponential relationship between them. The large porosity corresponds to smaller breakout pressure according to the fitting formula. In order to show the relationship between porosity (P) and permeability (p), their fitting profile is shown in Figure 7b. The fitting formula is P = 1.13 × 10 7 × e −8.71 p , which is a typical exponential relationship. In the process of oil and gas exploitation, the character of the pore structure is first acquired, the gas injection breakout pressure is calculated next by the intersection profile, and last, the permeability is predicted by log-log analysis in theoretical calculations.

Conclusions
The penetration and displacement behavior of molecules in interlayers with different porosities were analyzed using molecular dynamics simulations. The results show that interlayers with different pore sizes and porosities can strongly affect the distribution, diffusion, and displacement processes. The MSD and gas distribution show that N 2 can pass through the pores of interlayer structures and displace the H 2 O, Na + , and Cl − in the pore space into the vacuum phase. This shows that high porosity can result in highefficiency diffusion of N 2 , and high porosity in porous interlayers benefit the displacement of H 2 O, Na + , and Cl − . Interaction analysis between N 2 and H 2 O confirms that a strong interaction between N 2 and water can overcome the interaction between water and the pore wall to render rapid displacement of H 2 O, Na + , and Cl − . A log-log analysis and an intersection profile were presented to theoretically predict the relationship between permeability, porosity, and breakthrough pressure.

Data Availability Statement:
No new data were created or analyzed in this study. Data sharing is not applicable to this article.